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ABSTRACT 

The site survey for the Advanced Technology Solar Telescope concluded re- 
cently after more than two years of data gathering and analysis. Six locations, 
including lake, island and continental sites, were thoroughly probed for image 
quality and sky brightness. The present paper describes the analysis method- 
ology employed to determine the height stratification of the atmospheric turbu- 
lence. This information is crucial because day-time seeing is often very different 
between the actual telescope aperture (~30 m) and the ground. Two indepen- 
dent inversion codes have been developed to analyze simultaneously data from 
a scintillometer array and a solar differential image monitor. We show here the 
results of applying them to a sample subset of data from May 2003, which was 
used for testing. Both codes retrieve a similar seeing stratification through the 
height range of interest. A quantitative comparison between our analysis pro- 
cedure and actual in situ measurements confirms the validity of the inversions. 
The sample data presented in this paper reveal a qualitatively different behavior 
for the lake sites (dominated by high-altitude seeing) and the rest (dominated by 
near-ground turbulence) . 

Subject headings: atmospheric effects — methods: data analysis — site testing 
— telescopes 



1. Introduction 

The Advanced Technology Solar Telescope (ATST, Keil et al. 2003), with its 4-meter 
aperture, will be the largest solar telescope in the world. A suite of advanced instrumentation 
will allow for unprecedented resolution and precision imaging, spectroscopy and polarimetry 
at visible and infrared wavelengths. In order to address the most challenging problems in 
modern solar physics, ATST's highest priorities are oriented towards high-resolution and 
coronal science. Obviously, these are aspects in which the site selection will have a strong 
impact and therefore should be considered very carefully in the decision of the telescope 
location. The ATST Site Survey Working Group (SSWG) was created in 2000 with an open 



^The National Center for Atmospheric Research (NCAR) is sponsored by the National Science Founda- 
tion. 
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call to the international community. Its goal was essentially to test the best solar observa- 
tories in the world and produce a recommendation to maximize the scientific capabilities of 
the ATST. More especifically, a set of quantitative goals were set for the ATST site: 

• Clear time of 70%, or 3000 annual hours of sunshine. 

• 1800 annual hours of seeing with a Pried parameter ro better than 7 cm, including at 
least 100 continuous 2-hour blocks. 

• 200 annual hours of ro > 12 cm, including at least 10 continuous 2-hour blocks. 

• Large isoplanatic angle for adequate performance of the adaptive optics system over a 
large field of view. This goal is equivalent to demanding minimal turbulence at high 
altitude. 

• 480 annual hours with a sky brightness smaller than 25 miUionths of the disk value at 
a distance of 1.1 solar radii. The radial profile curve should follow a power law with an 
exponent i? = 0.8 or steeper. These conditions should be met in at least 40 continuous 
4-hour block per year. 

• 600 annual hours with precipitable water vapor below 5 mm, including at least 40 
continuous 4-hour blocks. 

Only two comparative studies of solar observatories had been conducted prior to this 
work. The JOSO organization studied nearly 40 sites in Europe to conclude that Observa- 
torio del Teide in Tenerife (Canary Islands) provided the better conditions recorded in their 
data (Brandt & Woehl 1982; Brandt & Righini 1985). The other work is that of the Caltech 
survey, which considered 34 sites in the southern California area to locate what is now the 
Big Bear Solar Observatory (BBSO; Zirin & Mosher 1988). 

Initially, a list of 72 well-established observatories and some promising undeveloped sites 
was compiled. Since the ATST SSWG had only sufficient resources for extensive testing of up 
to six sites, a first decision had to be made in the absence of any actual measurements. This 
decision was made based on previous observing experience of the group members. A first 
round of discussion quickly resulted in the selection of BBSO, La Palma (Canary Islands) 
and Sacramento Peak (New Mexico), which have an outstanding record in producing high- 
quality solar observations. The SSWG also selected Hawaii, and in particular Haleakala (on 
the island of Maui) was preferred over Mauna Loa and Mauna Kea after taking feasibility 
issues into account. Diverse other considerations later resulted in the choice of Panguitch 
Lake (Utah) and San Pedro Martir (Mexico) for the final two slots. 
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The instrumentation deployed at the six candidate sites includes a scintillometer array 
(SHABAR) and a Solar Differential Image Motion Monitor (S-DIMM). Both instruments 
have been described in a recent paper by Beckers (2001). The SHABAR and S-DIMM are 
used to determine the atmospheric turbulence as a function of height. It must be noted that 
daytime seeing conditions are often strongly dependent on height above the ground. Thus, 
the ATST site survey needed to determine the observing conditions at the telescope aperture 
height, approximately 35 meters above the ground (the actual height is still to be decided). 
The site survey instrumentation was placed at a height of 8 meters (6-meter stand plus a 
2-meter telescope pier) , which is high enough to avoid the lower part of the surface boundary 
layer but not so high that wind jitter could become a significant problem. 

Data from the SHABAR and S-DIMM instruments have been gathered and analyzed 
until August 2004. In October 2004, the ATST SSWG presented a final recommendation 
and produced an exhaustive document with a detailed report of the site survey process. The 
SSWG report is publicly available for download from the ATST website (http: / / atst.nso.edu/) 
That document contains abundant technical documentation and the reader is referred to it 
for further details. The present series of papers is a formal publication of the most rele- 
vant results. This first paper describes the inversion methods used to infer the atmospheric 
turbulence and the Fried parameter from the SHABAR/S-DIMM data. The methods have 
been extensively tested by means of simulations, comparison between different algorithms 
and also with in situ measurements. One month of data (May 2003) was used for testing 
and fine tuning the algorithms before being applied systematically to the entire site survey 
period. 



The final objective is to determine the Fried parameter as a function of height, ro[h). 
This parameter depends not only on the turbulence at a height h but also everything above 
it. With some assumptions (see Hickson & Lanzetta 2004 for a complete derivation), tq can 
be obtained as: 



where C is a constant, z is the zenith angle and (which represents the amplitude of 
the refraction index fiuctuations) characterizes the local atmospheric turbulence at height 
h. Notice that h denotes height above the instrument, which in our case is located 8 meters 
above the ground. 



2. Derivation of C„(/i) 



2.1. 



Basic relations 




(1) 
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The covariances {Bj) of the brigthness fluctuations measured by two detectors, separated 
by a distance d in the SHABAR array, are related to by: 



where the W{h, d) are known functions with analytical expressions. If a is the detector 
diameter and a is twice the tangent of the solar radius (which we consider approximately 
constant though the year): 



with Jo and Ji being the Bessel functions. The total scintillation can be obtained from 
the expressions above as Bi{d = 0). The kernels W{h,d) are shown in Fig 1. It can be 
deduced from this figure that the Bj corresponding to further separations (larger values of 
d) have their sensitivity peak at higher atmospheric layers. Fig 2 shows average values of 
Bi{d) measured at the six candidate sites. Not surprisingly, the covariances decrease with 
the detector separation d. 



The inversion strategy that we employed to solve for C^{h) in Eq (2) is a non- linear least- 
squares iterative fitting combined with singular- value decomposition (SVD) of the covariance 
matrix. The details of the implementation can be found, e.g. in Press et al. (1986), but the 
basic idea is the following. One starts with an initial guess for the C^{hi) (where hi are the 
gridpoints of the atmospheric discretization). With this initial model we calculate the 15 
different covariances Bi{djk) (where djk is the separation between each pair of detectors j and 
k in the SHABAR instrument), the total scintillation Bi{d = 0) and the Fried parameter at 
the base of the atmosphere ro{h = 0) (measured by the S-DIMM instrument), using Eqs (1) 
and (2). The derivatives of these parameters with respect to the Cf^{h) (the so-called response 
functions) are also computed, as this information is needed by the fitting procedure. The 
synthetic Bi and tq are then compared to those actually measured by the SHABAR and 
S-DIMM. In general, the synthetic values will not match the measurements. The algorithm 
then uses the information in the response functions to calculate the first-order corrections 
SC^{hi) that minimize the difference (in a least-squares sense) between the synthetic and 
measured values of Bi{d) and ro{h — 0). 




(2) 




(3) 



2.2. 



The inversion codes 
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Fig. 1. — Height dependence of the inversion kernels W{h,d) for the scintillation and the 
fifteen covariances (for various values of d) between brightness fluctuations measured by the 
SHABAR. The upper and lower curves correspond, respectively, to the smallest and largest 
sperations (0 and 468 mm), z is the solar zenith angle. 
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Fig. 2. — The average cross-correlation curves for the sites as of October 2003. The two 
lake sites, Big Bear (BB) and Panguitch (PG), have relatively shallow curves as a function 
of detector separation which indicates that the seeing arises from a region far above the 
instrument. The four non-lake sites have a relatively steeper cross-correlation as a function 
of separation, suggesting a boundary layer near the telescope. 
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The numerical solution of an inverse problem such as the one described here is always 
challenging, but even more is to ensure that the results obtained are physically meaningful. 
Understanding the issues and uncertainties involved in the inversion is critical for the success 
of the site survey and, ultimately, the ATST. It is equally important to ensure that the nu- 
merical codes employed are exhaustively tested and debugged. Based on previous experience 
with conceptually similar problems, the SSWG decided that two different codes should be 
independently developed and applied to the data. The two codes, based on the same min- 
imization scheme described above, exhibit some differences in the actual implementation. 
By comparing the results provided by two different codes, we can learn about the reliability 
of the procedure and the limitations of our seeing measurements. Since there is no reason 
to prefer one above the other, we take the difference between their respective results as an 
estimate of the error bar. 



The first code was developed at the Instituto de Astrofisica de Canarias (lAC) and run 
by National Solar Observatory (NSO) staff. The spatial grid is equispaced in log(/i) with 68 
points ranging from 20 cm to 40 km. However, not all of these points are inverted for (the 
problem would be ill-posed due to the large number of free parameters). Instead, C^ih) is 
determined at a subset of 17 special locations, the so-called inversion nodes. The rest of the 
atmosphere is then reconstructed via splines interpolation from the values at the nodes. 

In order to ensure that C^{h) is always positive, we used the following change of variable. 
The code solves for a variable y such that y{h) = log[C^^{h)]. In this manner, regardless 
of whether y is positive or negative, C^X^) will always be positive. Before inverting, the 
measured Bi and ro(/i = 0) are processed with a 5- minute running median filter to smooth 
out fluctuations due to instrumental issues or problems such as birds or insects flying by the 
detectors. 

After an initial application of the code to a sample dataset, it was found that a significant 
number of points could not be successfully fitted by the procedure. The observed ro{h = 0) 
provided by the S-DIMM instrument was inconsistent with the observed scintillation values 
Bj{d) from the SHABAR. In all cases the S-DIMM tq was smaller than would be expected 
from the scintillation measurements. To overcome this difficulty we introduced an additional 
parameter As in the inversion that accounts for this "missing scintillation", so that the 
synthetic Bj are now given by: 



2.2.1. The I AC code 




(4) 
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while Eq (1) is still used to compute ro{h). 

The necessity to include As points towards the existence of a a source of image degra- 
dation that does not produce scintillation. This may be due to a finite outer turbulence 
scale, to "windfiltering" (due to the finite time amount it takes for the turbulence to cross 
the field-of-view at high altitudes) or to high-altitude seeing (higher than ~1 km, where 

the kernels converge). In either case, the inclusion of As would be appropriate. While we 
have not been able to pinpoint the exact nature of the missing scintillation, the statistical 
distribution of As (see Fig 3) provides an interesting clue. 

The two island sites (La Palma and Haleakala) require a relatively small amount of 
missing scintillation, whereas the lake site (Big Bear) frequently demands a substantially 
larger value. As we discuss below, the seeing in the island sites is typically dominated by a 
boundary layer near the ground and the turbulence decreases rapidly as we move upwards. 
The lake sites, on the other hand, are frequently dominated by high-altitude seeing. Thus, 
it seems hkely that the As is related to the turbulence above the sensitivity range of the 



The second inversion code was developed independently at the High Altitude Observa- 
tory (HAO). One of the essential differences with the lAC code is that all the grid-points 
in the atmospheric discretization are now free parameters. In order to make the algorithm 
stable, a regularization scheme was implemented in such a way that the code prefers smooth 
solutions whenever possible. Basically, in each iteration we add a penalization term to the 
X^- This term is proportional to the quadratic deviation of C^{h) from a straight hue on 
a log-log plot. The fact that each grid-point is a free parameter represents an important 
burden that slows down the execution of the code. This is compensated, to some extent, by 
the more accurate integration scheme described below. 

The model atmosphere only goes up to a maximum height Hm = 1000/ sec{z) where 
the inversion kernels merge (Fig 1). We still consider the effect of turbulence above Hj^ by 
breaking down the integrals as: 



SHABAR. 



2.2.2. The HAO code 





and 
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Note that Bj does not depend on d. As before, a and Bj are free parameters that 
account for the scintillation and image degradation coming from high-altitude turbulence. 

Another important difference with the lAC code is the integration scheme used to solve 
Eqs (5) and (6). The integral over each finite grid interval is computed assuming a parabolic 
dependence of log(C^) with log(/i) (see the appendix for details). While this integration 
introduces some additional complexity in the algorithm, it also alleviates the computational 
expense of the inversions by allowing a coarser grid without loss of accuracy. 

The data pre-processing includes a 5- minute block average and a test to look for clouds in 
the averaged interval. If the derivative of the intensity measured by the SHABAR detectors 
changes sign two or more times, the entire period is flagged as cloudy and the ro(/i) is set to 
zero at all heights. 



3. Verification tests 

3.1. Compcirison of two codes 

A sample subset of data (those taken in May 2003) were subject to inversion with 
both codes. Fig 4 shows the cumulative distributions produced by each code at La Palma, 
Haleakala and Big Bear. The results are fairly similar for most of the height range represented 
in the figure. As one would expect, the discrepancy increases with height due to the limited 
sensitivity of the SHABAR measurements to turbulence at higher layers. The same behavior 
is observed with the in situ measurements (see section 3.2 below). 

A more detailed representation is given in Figs 5 and 6, which show the average strat- 
ifications of and Tq with height. The produced by the HAO code has more spatial 
structure (small-scale fluctuations) than the lAC code. However, the overall dependence is 
actually very similar and the lAC results resemble a smoothed-out version of those from 
HAO (at least down to the first meter above the detectors). Whether the fine structure 
inferred by the HAO code is real or not is practically irrelevant for our purposes here, as the 
ro(/i) parameter, which is our ultimate goal, is barely sensitive to it. The Tq stratifications 
from both codes arc very similar, at least up to ~ 40 m above the instrument. In this small 
sample we can already see the different behavior of the lake sites (Big Bear and Panguitch) 
as opposed to the island (La Palma and Haleakala) and continental sites (San Pedro Martir 
and Sac Peak). The Big Bear plot exhibits a larger tq near the ground, but it does not 
improve significantly with height. The island and continental sites, on the other hand, have 
a lower value at the ground but the seeing improves dramatically as we go to higher layers, 
suggesting the presence of a boundary layer near the ground. The slope is particularly steep 
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in the ~30-40 m range, precisely where the ATST aperture will be. These results emphasize 
the importance of analyzing the height dependence of tq, as we have done in this work, since 
the site and height of the telescope may end up being critical factors in the performance of 
the ATST. 



3.2. In situ measurements 

Perhaps the most rehable test of the ATST site survey analysis is to compare its results 
directly with actual in situ measurements. The SSWG had access to a tower with meteo- 
rological instrumentation including hygrometers and sonic anemometers at various heights 
(DASH, Oncley & Horst 2004). With this instrumentation it is possible to make local mea- 
surements of temperature T and humidity g up to a height of ~100 m. These parameters 
are monitored at a rate of 30 Hz, from which can be derived. 

A complete SHABAR/S-DIMM system was deployed in Erie (Colorado), at the base of 
the measuring tower. Simultaneous measurements from the tower and the ATST instrumen- 
tation have been used to produce the scatter plots in Fig 7. Each panel shows the diagonal of 
the plot, as well as the Pearson correlation between the ATST inversions and the SHABAR 
instrumentation. Before analyzing the discrepancies between the two, a brief explanation is 
in order. While the tower measurements are strictly local, the SHABAR detectors receive 
solar hght that has passed through a cone of Earth atmosphere. Thus, the ATST instru- 
mentation is sensitive to a field of view that becomes larger for higher heights. For this 
reason, the measured in situ is expected to exhibit larger temporal fluctuations than that 
from the ATST instrumentation, which is a spatial average over the SHABAR field of view. 
This effect, which can be seen in Fig 8, increases the scatter especially in the higher layers. 
Moreover, the SHABAR field of view is displaced with respect to the tower instrumentation 
by an amount that increases with height. In any case. Fig 8 shows that the overall trends 
are in general well retrieved by the ATST analysis. 

The tests indicate good agreement between both types of measurements up to a height 
between 22 and 50 m. It is important to note that the height range where the ATST 
measurements are accurate is lower in this case than at the actual candidate sites, simply 
because the tests were carried out during mid-winter at the considerably higher latitude of 
Erie (the sec [2;] angle was always larger than 2.0). 

While the Erie tests reassured our confidence in the ATST seeing measurements, there 
was still some concern that the atmospheric structure above the continental Colorado is 
probably very different from the actual sites. The SSWG then decided to mount similar 
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in-situ probes on a portable crane to perform similar tests on Big Bear and Haleakala. New 
measurements (this time considering lower heights due to the limitations of the crane) were 
conducted at the two sites, with the results shown in Fig 9. 

The Haleakala tests show excellent agreement between the inversions and the measure- 
ments. In Big Bear, however, the correlation is very weak. Among the differences between 
the two sites are the level of humidity, the local topography and environment, and the 
probable lack of a ground boundary layer in Big Bear. 

A possible explanation may be that the weaker scintillation in Big Bear makes these 
measurements more sensitive to noise. A closer look at the data shows that the correlation 
is improved by removing outlier points and clouds (Fig 10). Particularly interesting is the 
dataset gathered on May 12 2004, which reveals the dependence of the data with the wind 
direction (Fig 11). When the wind blows from the lake (west), the observing conditions are 
optimal and the seeing is dominated by high-altitude turbulence. On that particular day 
the wind changed from west to variable directions, resulting in the development of a ground 
turbulence layer and improving the correlation coefficient to ~0.77. Note also that, in the 
absence of a surface turbulence layer at Big Bear, tq will essentially be independent of h near 
the ground and the S-DIMM measurements at 8 m can be taken to represent ro at 25-40 m. 



4. Conclusions 

The ATST project has imposed stringent requirements on the telescope site, both in 
terms of seeing and sky brightness. The seeing analysis has proven very challenging due to 
the strong dependence of Tq with height exhibited by most sites. This paper describes the 
analysis employed by the ATST SSWG to find the site that provides optimal conditions (ro) 
at the height of the ATST primary aperture. 

Various tests have been used to demonstrate the validity of our analysis, at least up 
to ~40 m above the ground. The use of two independent inversion codes are helpful to 
understand the uncertainties in the procedure and to set error bars on the measurements. 
Moreover, in situ measurements proved the validity of the approach in the height range of 
interest. The only caveat identified by our work is that the errors increase when the seeing 
is dominated by high-altitude turbulence (i.e., in the absence of a ground boundary layer). 
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A. Pcirabolic integration scheme 

Let X = \og{h) and y = log[C^{h)]. If we approximate y and W as parabolic functions 

of x: 

W — WiX^ + W2X + W2 (Al) 

y — ax^ + bx + c, 

we can solve analytically the integral in Eq (5) between two adjacent grid-points Xi and X2, 
as: 



where, 



r [wix^ + W2X + W3] dx ( A2) 

10^ f ' [wix^ + W2X + ws] 10("^'+''^)dx 



= 10" (Ti + T2 + T3) 
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T, = ; /" ^ {2^+1^ 5^ ((10^1 (^+"^1) - 10^^ (^+«^^)) b 

54^ a2 log(10)2 

+2 a (- (I0^i(^+"^i)xi) +10^2{6+ax2)^^^^ ^log(lO) 

+V^Erf( ^^ + ^"";^/^ ) (2a - 6^ log(lO) ) 
2 v a 

+AErf( <'' + '''";'y^ ) (-2a + 6^ log(lO)) 

2 V a 
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54^ a2 log(lO) ^ 

2 ya 
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Fig. 3. — Histograms of the relative frequency of occurrence of missing scintillation As values 
at the three sites. Left: relative occurrence of As as a fraction of the observed scintillation; 
left: entire distribution; Right: relative occurrence of As as an absolute scintillation measure. 
Solid line: Big Bear; Dashed line: Haleakala; Dotted line: La Palma. Note that the two 
mountain sites (La Palma and Haleakala) typically require 10-20% fractional As while the 
lake site (Big Bear) frequently needs a substantially larger value. 
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Fig. 4. — A comparison between the two inversion methods. Here we show the cumulative 
distribution of the esimated rf){h) at the three sites and at four heights, as derived from the 
HAO method (dashed hne) and the lAC method (sohd hne). 
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Fig. 5. — Curves of C^{h) and ro(/i) retrieved by the inversions averaged over May 2003 for 
Big Bear, Haleakala and La Palma. Dashed: HAO method, Sohd: lAC method. 



- 19 - 



■3 



fit SPttIi 



L 



•.<ri 





T 

i 



t 1 



^3' 




-I i 4 g a 

I-ehiclIj:- UH* ryoA 

4ij • n' w H* ■ ■ ■ ■■■■■■ 4J ill 1 1 iBJ ■ ■ 



■ 




I- 
E 



vN" 



I' 



-1 



■■ ■■■ ■■■■ iiji ■■■■ 

E q 1 



II ■■ ill ■■ ^ 



■■■ ■ ■! ■■■ || I ^^^B- 




III ■■■ til III !■■ ■isi 

■ . „ J 

/ 

Ji 



1^ ■■■ ■■■ ■ ■■ ■■■ ■■■ 

c t t 



Fig. 6. — Curves of C^{h) and ro(/i) retrieved by the inversions averaged over May 2003 for 
Sac Peak, San Pedro Martir and Panguitch. Dashed: HAO method, Sohd: lAC method. 
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Fig. 7. — Comparison of from the ATST system with in situ measurements in Erie, 
Colorado. 
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Fig. 8. — Time series determinations of with the Erie tower and the ATST instrumenta- 
tion at 5, 10, 22, 50 and 100 m. Bottom-right panel: Median value as a function of height. 
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Fig. 9. — Scatter plot comparisons between the ATST estimates and the in situ measure- 
ments of at Haleakala (left four panels) and Big Bear (right four panels). In each four- 
panel set the comparison is for the three heights of the in situ measurements, plus all of 
the points combined. The straight line is strict equality, and each plot is labeled with the 
Pearson correlation coefficient of the data. 
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Fig. 10. — The influence of noise on tlie correlation between the ATST and in situ measure- 




Fig. 11. — The correlations for Big Bear on May 12, 2004 when the winds were atypically 
variable. Top: time series for the three heights. Bottom: correlation for this day. 



